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ABSTRACT 

We study the impact of relatively weak AGN feedback on the interstellar medium 
of intermediate and massive elliptical galaxies. We find that the AGN activity, while 
globally heating the ISM, naturally stimulates some degree of hot gas cooling on 
scales of several kpc. This process generates the persistent presence of a cold ISM 
phase, with mass ranging between 10^ and ^ 5 x 10^ Mq, where the latter value is 
appropriate for group centered, massive galaxies. Widespread cooling occurs where the 
ratio of cooling to free-fall time before the activation of the AGN feedback satisfies 
Aooi/te ^70, that is we find a less restrictive threshold than commonly quoted in the 
literature. This process helps explaining the body of observations of cold gas (both 
ionized and neutral/molecular) in Ellipticals and, perhaps, the residual star formation 
detected in many early-type galaxies. The amount and distribution of the off-center 
cold gas vary irregularly with time. 

The cold ISM velocity field is irregular, initially sharing the (outflowing) turbulent 
hot gas motion. Typical velocity dispersions of the cold gas lie in the range 100 — 
200 km s“^. Freshly generated cold gas often forms a cold outflow and can appear 
kinematically misaligned with respect to the stars. 

We also follow the dust evolution in the hot and cold gas. We And that the 
internally generated cold ISM has a very low dust content, with representative values 
of the dust-to-gas ratio of 10“^ — 10“^. Therefore, this cold gas can escape detection 
in the traditional dust-absorption maps. 

Key words: 

galaxies: elliptical and lenticular, cD; galaxies: evolution; galaxies: ISM; galaxies: 
groups: general. 


1 INTRODUCTION 

Elliptical galaxies generally have a complex multiphase in¬ 
terstellar medium (ISM), which is a key component of the 
galactic ecosystem. The mixture of the various gas phases is 
usually different in early-type galaxies (ETGs) with respect 
to late-type systems and this, in turn, reflects on the differ¬ 
ent evolutionary pattern of morphological diverse galaxies. 
Reviews of hot, warm and cold gas in ETGs can be found 


in e.g. Caon et al.j (2000 

); Mathews & Brighenti (2003a); 

Sarzi et al. (2006 

120131; 

Mulchaey & Jeltema (2010); El- 

lis & O’Sullivan ( 

2006); Young et al. (2011); 

Davis et al. 

(2011); Serra et al. (2012). The cold ISM is often spatially 


extended, with irregular distribution and kinematics (Caon 
et al.|2000 ). Dust is also present in all the phases of the ISM 
( Temi et al.|20Q4| 2Q09| Smith et al.|2012 ), and the dust-to- 
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gas ratio can be used to track down the origin of the cold 
gas (see below). 

A multiphase ISM is expected to keep the galaxy 
“alive”, allowing for (perhaps tiny) episodes of star forma¬ 
tion (SF). Indeed, a sizeable fraction of ETGs shows evi¬ 
dence of recent SF (Trager et al. |2Q0Q[ [Kaviraj et^20Q7| ). 
Therefore, a grasp on the origin and evolution of the multi¬ 
phase ISM would reveal vital information on the evolution 
of ETGs, such as star formation history, quenching mecha- 
nism(s), the role of active galactic nuclei (AGN) feedback, 
black hole growth, accretion of baryons from the environ¬ 
ment and more. 

In most ETGs a significant fraction of cold gas is 
thought to have an external origin, as indicated by the mis¬ 
aligned kinematics of the cold gas (both ionized and neutral) 
with respect to the stars (e.g. [Davis et al.||201l] and refer¬ 
ences therein). Cold gas is indeed commonly found around 
and bound to ETGs (Thom et al. 2012). The fate of this 


cold gas reservoir is still unclear, but it is likely that some 


© 2002 RAS 
























































2 M. Valentini & F. Brighenti 


of this low entropy gas is falling onto the galaxy and reaches 
the core in 1 Gyr. However, ACDM hydrodynamical sim¬ 
ulations on the kinematics of accreted HI in ETGs indicate 
significant disagreement with respect to the observations, 
making the source of misaligned cold gas still an open ques¬ 
tion ( Serra et al.|2014 ). 

Some of the cold ISM certainly comes from internal pro- 


demonstrate the physical existence and the basic properties 
of widespread gas cooling in ETGs. 

The thermal stability of gas in galactic or cluster cool¬ 
ing flows has been investigated since the discovery of the hot 
ISM and IGM, through both linear stability analysis and nu¬ 


merical simulations (see, for instance [Mathews Bregman 


cesses, such as stellar mass loss or hot gas cooling (Davis 
et al.|20rT Werner et al.[|20i4 David et al.|2014 ). Theoret¬ 


ical investigations about the formation of a cold ISM phase 


in ETGs include 

Brighenti & Mathews ( 

2002); Mathews & 

Brighenti (2003b 

); Temi et al. (2007); Parriott & Bregman 

(2008); 

McGourt et al. (2012); Gaspari et al.[([2012a|b ); Li & 

Bryan ( 

2014); Brighenti, Mathews & Temi (2015). Eor mas- 


msl 

Gowie et al.[[1980| [Bodo et al.|[1987 

[Malagoli et al. 

1987a 

Balbus 

1988a 

Loewenstein[ 1989a[ 

Balbus & Soker 

1989a 

Malago 

i et a , 

. [1990a[ [Hattori & Habe[[1990[ [Reale 

et al.[[ 

1991[ [Binney et al.[[2009[ [McGourt et al.[[2012[ [Joung 


sive galaxies, especially when located in high density envi¬ 
ronments, the internally produced cold gas may well domi¬ 
nate over the accreted one. Gavagnolo et al. (2008) and Voit 


et al.|2012 [Gaspari et al.|2012b Li & Bryan 2014). However, 
linear stability analysis provides an incomplete understand¬ 
ing when strong density perturbations are present. It is quite 
obvious that an overdense clump of gas can cool fast and sep¬ 
arate from the hot phase. The density contrast threshold for 
this localized cooling depends on the properties of the back¬ 
ground medium and the perturbation (see Joung et al.|2012 


et al. (2008) found that central galaxies in clusters host a 


multiphase gas only when the central IGM (ISM) entropy is 
below a threshold kT/n^^^ ~ 30 keV cm^, providing a clear 
link between the presence of cold gas and the properties 
of the surrounding hot gas. In a recent theoretical investi¬ 
gation Lagos et al. (2014) argue that cooling from twisted 
hot gas halos results in a large fraction (~ 46 %) of ETGs 
with internally generate d misaligned cold g as, explaining the 
ATLAS^^ observations ( Davis et al.|2011 ). 


In this paper we study the generation of a spatially 
extended cold gas phase in ETGs by radiative cooling of 
the hot halo gas. The key question is therefore if, and un¬ 
der which circumstances, the hot ISM can cool at large 
distance from the galaxy center. Multiwavelength observa¬ 
tions often show a strong spatial correlation between X-ray 


for a recent numerical investigation). Therefore, the key as- 
trophysical problem is to understand if such strong pertur¬ 
bations are naturally generated inside elliptical galaxies. As 
we discuss below, the answer delivered by our simulations 
is yes. We prefer to keep the assumptions of our models to 
a minimum. We do not use an ad hoc distributed heating 
to ensure global thermal equilibrium, nor we superimpose 
an artificial turbulent velocity field to generate density per¬ 
turbations. Adopting low-power, collimated outflows as the 
only perturbing means of the ISM, we aim to probe the im¬ 
portance of AGN-stimulated cooling in ETGs in a simple 
yet solid way. 


2 NUMERICAL PROCEDURE 


and Ha images (e.g. [Trinchieri et al.[[1997| Goudfrooij & 

The simulations have been carried out with a modified ver- 

Trinchieri[[1998[ [Trinchieri & Goudfrooij [2002[ [Temi et al. 

sion of the eulerian code ZEUS-2D (Stone & Norman|1992 ), 

2007 

Gastaldello et al.[[2009 Werner et al.[|2014), implying 

developed to include in the set of hydrodynamical equations 


rapid local cooling. Eurthermore, the connection between 
cold gas, radio emission and AGN-driven X-ray disturbances 
(see references above, also Blanton et al.|20rT ) suggests the 


source and sink terms describing galactic scale cooling flows. 
These terms include mass and energy injection by the stellar 


key point that spatially extended cooling has been triggered 
by the AGN activity. Evidently, while AGN feedback sup¬ 
presses the total cooling rate with respect to standard cool¬ 
ing flow models (e.g. McNamara & Nulsen [2007 ), it also 


population, as described in Mathews & Brighenti (2003a) 


stimulates some gas condensation on scales of several kpc. 

Guided by these observational results, we present here 
a detailed investigation of spatially extended gas cooling re¬ 
sulting from AGN feedback. We do not address the long term 
fate of the cooled gas — it can be in part accreted onto the 


o, T X • (p d u) — QxcP* (5* Psputt H“ Pgrowth 5 (l) 


central black hole (Gaspari et al. 2012b), it can settle in a ro¬ 
tating disk ( Brighenti Mathews|19% 1997a) or it may be 
re-heated to the hot gas temperature by thermal conduction 
or other processes. 


The present time type la supernova rate adopted here is 0.02 
SNrQ The code has been further enhanced to deal with a 
second fluid representing the dust. Dust density pd obeys to 
the equation: 

d Pd 
dt 

where u is the velocity of the gas, a* the specific rate of 
stellar mass loss, p* the stellar density and a uniform stellar 
dust to gas mass ratio (5* = 0.01 has been assumed. On the 
right hand side of equationbesides the stellar source term 
for dust, there are the sink term due to grain sputtering and 
the source term due to grain growth, psputt and Pgrowth, 
respectively. 


Eor the reasons outlined in Gaspari et al. (2012a) we 
consider collimated, non-relativistic outflows as the domi¬ 
nant feedback mechanism in local ETGs. The AGN feedback 
process is far more complicated than we consider here and is 
not appropriately understood (see, e.g. [Gaspari et al.|2012a 


Thermal sputtering (Draine & Salpeter 1979 Tsai & 
Mathew^[1995 ) rate is estimated as: 

pd 


Psputt 


^sputt 


= rid (mgr) , 


( 2 ) 


for a discussion on the numerical limitations in simulating 
AGN heating in galaxies and clusters). However, for the ob¬ 
jectives of this paper it is not necessary to accurately model 
the whole feedback cycle. An idealized feedback scheme as 
the one described in the next section is adequate enough to 


where rid — Pd/(^gr) is the number density of dust grains. 


^ These terms have no significant impact on our results and have 
been included just for consistency with previous cooling flow cal¬ 
culations. 
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(rhgr) is the average mass sputtering rate per grain, rhgr = 
A'ko? pgr a being the mass sputtering rate for a grain of given 
radius a. Here, pgr = 3.3 g cm“^ is the density of silicate 


grains ( Temi et al.|2003 ), a is the rate at which the grain ra¬ 


dius decreases ( Tsai &; Mathews|1995 Mathews & Brighenti 
2003b) and the resulting average sputtering time is: 


"T’sputt — 


(rhgr) 3|d| a‘^~^ da 3|d|(a2) 


0.03 


1 + 


2 X 10^ 


6 \ 2 . 5 ' 


(3) 


Myr . 


In equation]^ s = 3.5 ( [Mathis et al.|1977| , amin = 0.01 /xm 
and Gmax — 1 Finally, (a^) = f a^~^da/ f a~^da and 

similarly (a^). 

The grain growth rate in cold and dense gas (jDwek 
1988| [Hirashita &: Kuo||2011| |Hirashit^|2012| |Zhukovs~ 
et al.||2008|) is estimated as: 


Pgrowth — ric^ (rhgr, growth) 5 


(4) 


where the growth rate of each grain of mass m and radius a is 
rhgr, growth = TTO^ fnmrrimVm, f = 0.3 being the probability 
that a colliding atom/metal sticks on the grain surface, rim 
and Vm oc a/T ( Dwek|[l998 ) the metal number density and 
thermal velocity (assuming an atomic weight mmlrrip = 20 ). 
The averaged dust density growth rate is therefore: 

A _(rhgr, growth) _ 3 ppd (« ) 

Pgrowth-/5d 

where p is the gas density and it has been taken into account 
that when the dust to gas mass ratio 5 approaches the value 
0.01 all the relevant metals are locked in dust grains. 

Radiative cooling describing the loss of energy for X- 
ray emission has been included with the term —nirreA(T) 
in the thermal energy equation, where the cooling function 
A is calculated according to Sutherland k, Dopitaj ( |1993| ), 
assuming solar metallicity. The cooling has been extended 
to tempera t ures l ower than T = 10^ K following Dalgarno 


& McCray (|1972| and assuming a value x = 10 of the 


fractional ionization. The cooling is truncated at T = 50 K. 
Dust-induced cooling, important for T > 10® K for a dust- 
to-gas ratio ~ 10 “^, has been taken into account according 
to Mathews & Brighenti (2003b). We are well aware that the 


cooling and chemistry of cold gas is extremely complicated 
and our treatment for the T < 10^ K gas is approximate. 

The two-dimensional z x R computational domain is 
made up of 1100 x 700 zones. The resolution is Az = AR = 5 
pc up to ^ X R = 5 X 3 kpc. The outer 100 points in both 
dimensions are separated by increasing distances which are 
the terms of a geometric progression whose common ratio 
is 1.095. This peculiarity of the grid allows to study the 
central region of the computational domain, where most of 
the cooling process is expected to occur, in an accurate way. 


bright group, the isolated galaxy {IG)^ an isolated, fairly 
massive elliptical, and the low mass galaxy (LM), an inter¬ 
mediate elliptical, which represents a more common member 
of the early-type galaxy family. The key difference between 
the first and the two latter models is the amount of hot 
gas within the galaxy, which in turn reflects on a significant 
difference in the cooling time profile (see Figure]^. In the 
CG most of the hot gas within 10 kpc is expected to be in¬ 
tragroup medium (IGM) (e.g. Brighenti Mathews||1999 ), 
while for the IG and LM all the hot gas originates from stel¬ 
lar mass los s ([Mathews k, Loewensteinj 19^ Loewenstein 


fc Mathews[[1987[ [Ciotti et al.[[1991| [Mathews fc Brighenti 

2003a[ ). 


2.1.1 Gentral Galaxy 

We model the GG to agree with NGC 5044, the bright¬ 


est galaxy in the homonymous group (see Gaspari et al. 


2011b 2012a for further details). The gravitating mass con¬ 
sists of a NFW ( Navarro et al. [P'996 ) halo with virial mass 
Mvir = 4 X 10^® Mq and concentration c = 8.5, plus a de 
Vaucouleurs prohle ( Mellier Mathe'^[1987 ) with effective 
radius and stellar mass Xe = 10 kpc and = 3.5 x 10^^ 
M 0 . Although inconsequential here, we also include a cen¬ 
tral black hole with mass Mbh = 4 x 10 ® M©. Regarding the 
hot gas, we adopt the observed temperature profile for NGG 
5044 ( Buote et al.[[2003l 2004 David et al.[[2009 ) as initial 
temperature distribution, then we calculate the gas density 
profile assuming hydrostatic equilibrium, choosing the cen¬ 
tral density to approximately agree with observations. Be¬ 
fore to activate the AGN outflow, we let the system devete 
a classical cooling flow, by evolving the flow for 100 Myi[j 
During this time, the initially present gas mixes with gas 
shed by the (old) stellar population of the galax>j^ The gas 
density and temperature profiles just before the activation 
of the AGN outflow are shown in Figure ^ The hot gas 
mass within r = 10 kpc is ~ 2.6 x 10 ^ M©, and the bolo- 
metric X-ray luminosity within the virial radius (rvir ~ 900 
kpc) is Lx ~ 1.5 X 10^® erg s“^. In Figure]^ it can be seen 
where this galaxy model is located in the B-band luminosity 
versus the bolometric X-ray luminosity plane for early type 
galaxies. As expected, it resides in the high luminosity tip 
of the observed Lb — Lx relation ( Ellis &: O’Sullivan|2006 ), 
populated by massive galaxies at the center of groups or 
sub-clusters. 


2.1.2 Isolated Galaxy 

For this model, representative of massive but non central 
ellipticals, the stellar profile follows a de Vaucouleurs dis¬ 
tribution, with total mass M* = 5 x 10^^ Mq and effective 
radius re = 6 kpc. The dark halo has a NFW profile with 
virial mass Mvir = 10^® Mq and concentration 10.0. 


2.1 Models 

We consider three models of elliptical galaxies: the eentral 
galaxy (thereafter GG), located at the center of an X-ray 

^ We checked that these values of minimum and maximum grain 
size do not affect our results. 


^ This has a negligible effect on the results presented in the fol¬ 
lowing. We decided to kick off the AGN feedback from a pure 
cooling flow just to be co nsistent with the IG an d LM m odels. 

^ We refer the reader to [Mathews fc Brighenti[ ( [2003a[ ) for the 
hydrodynamical equations solved here and the relevant source 
terms for mass and energy. 


© 2002 RAS, MNRAS 000, [ip^ 






































































4 


M. Valentini & F. Brighenti 
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Figure 2. Lb — Lx diagram for the elliptical galaxies sample 
of |Ellis fc 0’Sullivan| ( [200^ . Filled circle represent galaxies with 
detected X-ray emission; open circles are used for upper limits. 
The location of our models, just before the AGN activation, is 
marked with the red squares when only the X-ray luminosity of 
the hot gas is considered. Big blue dots individuate the models 
when the contribution of the low mass X-ray binaries is added 
( |Kim &: Fabbiano|20Q4| >. 


Figure 1. Upper panel: Density profile just before the outflow 
activation for the central galaxy (solid line), the isolated galaxy 
(dashed line) and the low mass galaxy (dotted line). The points 
represent data for NGC 5044 taken from |Buote et ah'] ( |2003[[20041 ); 
[David et al.| ( [20Q^ . Bottom panel: initial (3D) temperature pro¬ 
files for the central galaxy (solid line), the isolated galaxy (dashed 
line) and the low mass galaxy (dotted line). Data points for NGG 
5044 are taken from jBuote et al^ ( |2003|[2004| ); [David et al.| ( |2009| >. 


The initial ISM is generated by evolving a classic cool¬ 
ing flow for 1 Gyr (the choice of this time is arbitrary and 
has no consequence on the results, since the flow reaches 
a quasi-steady state after ~ 100 Myr). The ISM entirely 
builds from the mass loss from the stellar population (e.g. 


Loewenstein & Mathews 1987 Mathews & Brighenti 20Q3a[ ). 
Just before the AGN outburst, the hot ISM mass within 10 
kpc is ~ 1.9 X 10® M©, the X-ray bolometric luminosity 
being Lx ~ 1.2 x 10^^ erg s“^. Figure shows that this 
model occupies an intermediate and uninteresting place in 
the Lb — Lx diagram. The initial gas density and temper¬ 
ature profiles are shown in Figure as dashed lines. Notice 
as pure cooling flows in elliptical galaxies show an overall 
negative temperature gradient, because of the central steep 
potential well and the resulting compressional heating. 


2.1.3 Low Mass Galaxy 


This model is a scaled down version of the IG described 
above. The stellar mass is now M* = 1.5 x 10^^ M© and the 
effective radius — 3.16 kpc, a value which agrees with the 


magnitude-re relation ( [Faber et al.[[l997 ). The NFW dark 
halo has a virial mass Mvir = 1.5 xTo^ M© and concen¬ 
tration 11.9. The mass of hot gas within r < 10 kpc is only 
~ 3.4 X 10^ M© (~ 1.15 X 10^ M© within re), with a X-ray 
(bolometric) luminosity Lx ~ 10^° erg s“^, which places this 
model within the lower envelope of the Lb — Lx relation. 

As common in intermediate and low mass ellipticals, 
only the gas in the inner region is inflowing, while the ISM 
at larger radii forms a subsonic outflow, powered by the 
SNIa (e.g. Mathews & Baker 1971[ [Mathews fc Brighenti 
2003a). For the assumed SNIa rate (0.02 SNu), the transition 
between inflow and outflow is located at r ~ 2 kpc. We note 
that for this relatively small galaxy the properties of the ISM 
are quite sensitive to the SN heating. Doubling the SNIa rate 
results in a global outflow (but the inner ~ 200 pc), with 
velocity increasing almost linearly from ~ 0 at 200 pc to 
~ 60 km s“^ at 3 kpc. The hot gas mass (within 10 kpc) is 
reduced to 1.9 x 10^ M©, while the X-ray lumi nosity drops 
to a very small value of 1.8 x 10®^ erg s~^ (see Mathews & 
Brighenti (2003a) for more details on the effect of the SNIa 


rate on the gas flow in elliptical galaxies). 
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Figure 3. Scheme showing the method we chose to implement a 
conical outflow (see text). The inner corner of the 2D cylindrical 
grid, with the 100 X 100 pc^ active region highlighted, is repre¬ 
sented. We note that in the code the 2 ; and R components of the 
velocity are located on different positions (see|Stone & Norman| 
Here we show the velocity vectors originating from the grid 
center just for illustrative purposes. 


2.2 Modeling the outflow 


2.2.1 Single event 


A subrelativistic and massive outflow (Brighenti & Mathews 
Gaspari et al.|[2012a ) which lasts 2 Myr (e.g. Fabian 


2006 


2012) is our single event AGN mechanical outburst. The out¬ 


flow energy injection is simulated dispensing velocity to the 
gas located in a box of {Rx z) = (100 x 100) pc. This aetive 
region is the inner corner of the computational domain, near 
the origin of the R and z-axi^ 

The number of cells by which resolve the interior of 
the jet has been chosen arbitrarily, the reasonable assump¬ 
tion is that the physics of the outflow does not influence 
the AGN feedback effects on kpc scale-length and the IGM 
global properties. 

The outflow velocity is given to the cells of the active 
region while the AGN is switched on, this energy input rep¬ 
resenting a spatially localized source term in the hydrody- 
namical equation of momentum. The velocity is given setting 
its module and defining a way for describing the components 
along R and z axes. We simulate both cylindrical and con¬ 
ical outflows, although we will focus mainly on the latter 
case. Cylindrical outflows have been simulated by equating 
the z-component of the velocity to the velocity module and 
setting the R-component equal to zero in all the aforemen¬ 
tioned zones. On the other hand, for conical outflows we 
set the velocity vectors as illustrated in Figure At any 
given grid point the velocity is radially directed, with the 
center located at (z,R) = (0, —100 pc/tg(30°)). The results 
are quite insensitive to the exact method of implementing a 
conical outflow. 

All the models have been evolved for 50 Myr, in or¬ 
der to investigate the short term evolution of a single AGN 
outburst. Each performed simulation was given a name: the 
identifying convention encodes the model of galaxy and both 
the shape and the velocity in km/s of the outflow. For in- 


^ Numerical experiments show that injecting the kinetic energy 
from a nozzle along the inner boundary conditions at 2 ; = 0 does 
not change the results. 


stance, model IG 4 OOO refers to the conical outflow with 
Vjet = 4000 km s“^ activated in the isolated galaxy. 


2.2.2 AGN feedhaek: reeurrent outbursts 


In these models we consider the whole AGN feedback cycle, 
where the gas cools and accretes on the central black hole, 
triggering an AGN outburst which heats the ISM. When 
the gas cools in the central region, the cycle is assumed to 
start again. Needless to say, the physics of AGN feedback is 
very complex and poorly known. We have no pretension to 
provide a self-consistent model here, we just need a simple 
scheme to set off outflows intermittently, in order to mimic 
the dynamical — sometime violent — evolution of the ISM 
in elliptical galaxies. We thus employ a simplified version of 


the feedback scheme used elsewhere (e.g. Brighenti & Math- 
ews||2006 Gaspari et al.|[^I2a| . Although this feedback is 
not fully self-regulated, the modeled outflows reproduce rea¬ 
sonable energies and powers (see Section |^. 

The adopted recipe has two main ingredients: the way 
of numerically describing the outflow and its activation tim¬ 
ing. We adopted a simple method in which the mechanical 
feedback activation is triggered by the presence of cooling 
gas in the innermost region (r ^ 100 pc) of the simulated 
galaxy. This assumption mimics the so-called cold accretion 
idea (e.g. Gaspari et al.|20I2b| and references therein): AGN 
outflows are fed by accretion of gas cooled out of the hot 
ISM, which accretes onto the central black hole. 

After several numerical experiments, we decided 
to adopt a drop out term —q(T)p/tcooi, w ith q{T) = 
qo exp[—{T/Tc)^], ( Brighenti Mathews|2002 ) in the conti¬ 
nuity equation to individuate and remove from the grid gas 
currently cooling to low temperatures. When gas dropout 
within r ^ 100 pc occurs (at least I M© in one timestep), 
the outflow is set off. We choose go = 2 and Tc = 5 x 10^ 
K to describe g(T), although the results are insensitive to 
the exact values of these parameters. This dropout method 
might result unpalatable and contrived, but has been shown 
to provide realistic estimates for cooled gas masses, limiting 
the effect of numerical overcooling over the traditional ap¬ 
proach (see [Brighenti Mathews ( |2QQ2| , Brighenti, Math¬ 
ews & Temi (2015) and the Appendix). Beside rectifying the 
numerical overcooling problem, there is an even more crucial 
reason to prefer the dropout technique. In order to calculate 
the genuine cooling rate at large distances from the center, 
it is preferable to avoid the presence of cold gas in the in¬ 
ner region, which would be ejected through an AGN outflow 
possibly reaching the off-center region and mixing with the 
freshly cooled hot gas. 

We note that we use the dropout term only near the cen¬ 
ter of the system (for r < 500 pc), so that the off-center gas 
cooling (see below) is calculated in the standard way, that is 
letting the gas cooling to low temperatures and keeping it on 
the numerical grid, although this is know to cause spurious 
overcooling (see Appendix). 

Once the outburst has been activated, the gas which 
resides in the inner region defined by R < 100, 2 ; < 100 pc 
is given a constant velocity as long as the condition on the 
presence of dropping out gas in r ^ 100 pc is met. 

All the models consider a sequence of conical outflows, 
their opening angle being 0 = 30°, and evolve up to 50 Myr. 
Models of recurrent feedback (FB) were given a name which 
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shows the model of galaxy and the velocity in km/s of the 
outflows, for instance CG-FB2000. 


3 COOLING AND TIMESCALES 

Here we briefly discuss the relevant timescales for thermal 
instability for our three models. It is well known that linear 
(infinitesimal) perturbations in the central regions of ETGs 
and galaxy clusters are generally stable in the pure hydro¬ 


dynamical case (e.g. [Malagoli et al.j 

|1987b| 

[Balbus [1988b[ 

Loewenstein||1989b| [Balbus & Sokerj 

|1989b| [Malagoli et al.[ 

1990b| Binney et al.||2009). In fact, slightly overdense fluid 


elements oscillates at the Brunt-Vasaila frequency ubv — 
27r/TBv, where usually the oscillation period tbv ^ Gooi- 
The oscillation center of the perturbation slowly sinks to¬ 
ward the center of the system, in a time only slightly shorter 
than the flow time of unperturbed gas. 

However, the linear stability analysis is of limited use, 
as in realistic galactic environments violent dynamical pro¬ 
cesses are expected to generate perturbations of finite ampli¬ 
tude. Indeed, X-ray images of cool cores clearly show large 
variations of gas density on all scales probed by current tele¬ 
scopes. 


Brighenti & Mathews (2002) have shown that spatially 


extended gas cooling in a cool core can happen when the 
hot gas is perturbed by the action of AGN feedback heat¬ 


ing. Brighenti & Mathews (2002) found that localized non¬ 


linear perturbations form in regions of converging flow and 
are able to cool without undergoing oscillations. Recently, 


McGourt et al. ( 2012 ); Sharma et al. ( 2012 ); Gaspari et al. 


( 2012 b); Singh & Sharma (2015) have studied in detail the 


cooling process in galaxy cluster cores perturbed by AGN 


feedback (see also Gaspari et al. 2011a). They provide a 
heuristic criterion for the onset of off-center thermal insta¬ 
bility: the IGM undergoes widespread cooling when stirred 
by AGN feedback if tcooi/^ff ^ 10, with the free fall time 
ts ~ tdyn = { 2 rlgY^‘^^ g being the local gravitational ac¬ 
celeration. [joung^^T^T ( 2012 ) investigated the evolution of 
nonlinear perturbations in the Galactic halo and found that 
for tcooi/Gcc ^ 1 the “cloud” cools before being disrupted 
by hydrodynamical instabilities. Here tacc = Cs/g, which is 
similar to tff. 

In Figure]^ we show the ratios tcooi/Gcc and tcooi/^ff for 
our three systems, at the time just before the AGN activa¬ 
tion. For the GG model the instability criterion tcooi/^ff ^ 10 
is met for r < 6 kpc, while for IG and LM no unstable re¬ 
gion exists outside the galactic core r < 100 pc, according 
to 


Sharma et al. ( 2012 ). 


However, in the following we show that the commonly 
adopted instability criterion is too restrictive. We find 
that extended cooling within elliptical galaxies occurs for 
Gooi/tff < 70 for realistic, recurrent, gentle AGN outbursts. 
We note that the empirical entropy criterion for star for¬ 
mation and multiphase gas in galaxy clusters cool cores. 


K = kT/n^/^ < 30 keV cm^ (Voit et al. 


2008 


Gavagnolo 


et al.|[ 200 ^ describes more accurately the spatial extent of 


the region where cooling happens. For the profiles shown in 
Figure we find that the radius where A = 30 keV cm^ 
is 17, 5.4, 3.6 kpc for model GG, IG and LM, respectively. 
These numbers compare well with the maximum size of the 
region where cooling occurs. 


4 SPATIALLY EXTENDED COOLING: SINGLE 
AGN OUTFLOW 

Below we analize the dynamics and the off-center cooling 
process triggered by a single AGN outburst. We focus here 
on the conical outflows, the cylindrical ones giving quali¬ 
tatively similar results. This idealized case will be used to 
better understand the physical nature of the non-linear den¬ 
sity perturbations which lead to widespread cooling. We ex¬ 
pect that in a more chaotic and realistic ISM, perturbed by 
the intermittent action of AGN feedback, extended cooling 
would be facilitated. This is the subject of SectionIn the 
following we focus on the gas that cools in the “off-center 
region” defined as: > 500 pc, R > 150 pc. We exclude the 
region near the symmetry z-axis because the 2D cylindrical 
grid used in our simulations might lead in spurious enhanced 
cooling along the axis. Indeed, the formation of cold gas 
along the jet axis is a common result of our simulations. This 
is physically reasonable because when the cavities generated 
by the AGN jets move buoyantly toward larger distances, 
they trigger a vortex flow around them, which compresses 
the hot gas toward the z-axis, significantly lowering the local 
cooling time (see also Mathews & Brighenti (2008) 


Revaz 


et al. ( |2008| ; Brighenti, Mathews & Temi(2015)). Prelimi¬ 


nary 3D calculations verify the formation of the z—axis cold 
filament (G. Melioli, in preparation), but we still prefer to 
be conservative in this work and focus on the region far from 
the symmetry axis for most of our analysis. 

We also neglect the gas which cools at the very cen¬ 
ter. The suppression of the total cooling rate in galaxies and 


clusters relates to the so-called cooling flow problem (Mc- 
Namara Nulsen|20Q7 ) and will not be discussed here. For 


calculations about AGN heating in ellipticals see [Gaspari] 
et al. (e.g. 2012 a), which use a similar feedback scheme. 


4.1 AGN-cooling in the Central Galaxy 

In the top panel, central column of Figurej^we show the time 
evolution of the cold gas mass located off-center {z > 500 pc, 
R > 150 pc), for the single conical outflow models. Times 
are measured since the activation of the outflow. For the 
GG, the three models with I’jet = 2000 , 4000, 8000 km 
s“^ are able to cool a substantial amount of gas, 10 ^ — 10 ^ 
M 0 , at large distances from the center. A weak trend with 
I’jet is evident from Figure]^ the larger is I’jet, the larger 
is the mass cooled off-center. As already found in [Brighentij 
& Mathews ( 2002 ), the AGN feedback generates non-linear, 
almost isobaric perturbations in regions of sustained com¬ 
pression, which trigger intermittent and widespread cooling 
throughout the region with tcooi/^ff ^ 10. Because of the 
large overdensities formed in converging flows, it seems im¬ 
proper to describe this cooling process as a thermal insta¬ 
bility. This agrees with recent studies of cooling in galaxy 
clusters by McGourt et al. ( 2012 ); Gaspari et al. (2012b). 


Gontrary to these latter works, however, we do not include 
an artificial source of turbulence nor a distributed source of 
heating. Generally, the region over which cooling occurs is 
more extended when the power of the outflow is larger. 

Interestingly, cooling can occur long after the AGN out¬ 
flow ceases (see also Brighenti, Mathews & Temi, 2015). For 
instance, for models GG 4 OOO and GG8000 off-center cooling 
episodes happen at t ~ 35 Myr and t ~ 50 Myr, respectively. 
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Figure 4. Upper panels: Cooling time, acceleration time and free fall time (see definitions in the text) just before the AGN outburst, 
for the CG (left), the IG (middle) and LM (right). Bottom panels: ratios of relevant timescales. 


while the outflow terminated at t = 2 Myr. At these times, 
usually no sign of the AGN feedback (cavities, shocks) would 
be detected with X-ray observations. Nevertheless, subsonic 
motion with velocity up to ~ 1/2 of the sound speed is still 
present in localized regions, especially in form of relatively 
large eddies with typical size ~ 0.5 kpc. We note that these 
features are well resolved in our simulations with ~ 50 — 100 
grid points. When exceptionally large eddies are present in 
the inner region of the galaxy, they appear as weak X-ray 
cavities. This happens infrequently in our simulations, with 
large outflow velocities and earlier times making this event 
more likely. 


The spatial distribution of the cold gas in the off-center 
region at several times, for models CG 4 OOO and CG8000, 
is displayed in Figure The cold blobs generally arise in 
a region close to the outflow symmetry axis {R < 1 kpc) 
and for 2 : < 10 kpc. The cold gas is arranged in a large 
number of small clouds (actually, toroidal structures, given 
the imposed 2D cylindrical symmetry), generally few 10s pc 
in size. As expected (see Koyama Inutsuka||2004 ) this is 
close to our numerical resolution (Ai? = /\z = 5 pc), with 
the unfortunate implication that the evolution of the cold 
clouds may be subject to serious numerical errors. We ad¬ 
dress some of the numerical problems, such as overcooling, 
in the Appendix. The off-center cooling process starts in a 
very localized fashion, usually in just one or few zones (see, 
for example, the top-left panel in Figure]^. This is true for 


every resolution adopted in our numerical experiments. Cal¬ 
culations using coarse resolution might be grossly in error. 

Using a tracer passively advected with the flow, we are 
able to track down the original location of the cooling gas. 
We find that most of the cold gas found off-center comes 
from the inner region (r < 1 kpc), where the entropy is lower. 
This confirms the results by |Li Bryan (2014). Because the 
ISM (iron) abundance peaks in the center, it also follows that 
the metallicity of the cold gas is somewhat higher than the 
one of the surrounding hot gas. 

In Figure we also show, in the central row, the evolu¬ 
tion of the cold gas mass in the region r = \^R‘^ > 500 

pc (that is, we include here the cold filament along the 
z—axis) and the total cold gas (bottom panels). The mass 
in the filament forming in the wake of the outflow lies in 
the range Mm ~ 10^ — 10® M(^ although these values must 
be confirmed by high resolution 3D simulations. The time 
evolution of the total cold gas mass shows that the single 
AGN outburst is unable to significantly suppress the cool¬ 
ing rate with respect to the pure coling flow scenario. Only 
for the GG8000 model the AGN outflow temporarily reduces 
the cooling rate by a factor of ~ 2.5 for a time ~ 15 Myr, 
after which the system recovers the standard cooling rate, 
M ~ 10 Mq yr“^. 


® Note that here and in the following we consider only the 2 ; > 0 
volume. 
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Figure 5. Upper left panel: Time evolution of the cold gas mass in the off-center region (R > 150 pc, 2 ; > 500 pc), for CG models 
with cylindrical jets CGeyl2000^ GGeyl4000 and GGeyWOOO. Here and in the following, times are measured since the start of the AGN 
aoutflow. Middle left panel: as above for the total cold gas mass in the region r > 500 pc, which includes the filament along the 2 ;—axis. 
Lower left panel: same for the total cold gas mass, most of which is located at the very center of the galaxy (r < 300 pc). Gentral panels: 
same as left panels but for the GG galaxy with conical jets (30 degrees) models GG 2000 , GG 4 OOO and GG8000. Right panels: same as 
central panels but for the isolated galaxy. 


4 . 1.1 Cylindrical outflows 

We have also calculated several models with cylindrical out¬ 
flows, in which the velocity is directed along the z—axis. The 
time evolution of the cold gas mass for the models with jet 
velocity 2000, 4000 and 8000 km s“^ is shown in the first col¬ 
umn of Figure The results are qualitatively similar to the 
conical outflows, with a tendency to generate more off-center 
cold gas. The reason is that cylindrical outflows are less effi¬ 
cient to distribute their kinetic energy in the inner region of 


the galaxy, therefore a larger volume of the system hosts an 
almost unperturbed cooling flow. Also, it is well known that 
conical AGN outflows get collimated by the hot gas pressure 
as they propagate outward, becoming effectively cylindrical 
at some distance from the center (cfr. [Brighenti Math- 


ews||2006 ). Given the similarity of models with conical and 


cylindrical outflows (and the uncertainties about real AGN 
outflows), in the following we limit the discussion to conical 
jets. 
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Figure 6. Left column: temperature map for model CG 4 OOO at three different times: t = 14, 32, 36 Myr, from top to bottom. The z-axis 
is horizontal, the i?-axis is vertical. Units are kpc. Time is measured since the switch on of the AGN outflow (which lasts 2 Myr). Right 
column: same for model CG8000. The temperature is shown at t = 18, 20, 22 Myr, to follow a strong cooling episode. Notice that the 
color scale differs for every panel. 


4.2 AGN cooling in the Isolated Galaxy 


As pointed out in Section the fundamental difference be¬ 
tween the IG and the CG models is that in the former the 
whole hot ISM is generated by the stellar mass loss of the 
galactic (old) stellar population. Therefore, the mass and 
the average density of the hot gas is much lower than in the 
GG case. As a result, the ratio tcooi/^acc for the IG is every¬ 
where larger by a factor of > 10 than for the GG (see Figure 
4|. According to the criterion by [Sharma et al.| ( |2Q12| ) and 
McCourt et al. (2012) we should not expect spatially ex¬ 
tended cooling for this model. This conjecture does not pass 
the test of our simulations. In fact the model IG2000 gives 
rise a substantial mass of cold gas at t ~ 5 — 25 Myr (top- 
right panel of Figure]^. We will see below that in case of 
repeated AGN outbursts, even moderate outflows with ve¬ 
locity of 1000 km s“^ trigger significant off-center cooling. 
The spatial map of the cold gas for IG2000 is shown in Fig¬ 
ure The off-center cooling process initiates at t ~ 5 Myr, 
with small cold blobs forming at {R,z) ~ (0.2,2) kpc and 
{R,z) ~ (0.4, 2.1) kpc, where the time ratio tcooi/^acc ~ 60, 


just before the AGN activatiorQ As time goes by, more gas 
cools and moves away from the original formation location. 
At later times, as shown in Figure a multitude of cold 
clumps can be found for {R, z) < (0.7, 2.5) kpc. As for the 
GG simulation, the emergence of the off-center cooling de¬ 
pends on two effects. First, the outflow transports (low en¬ 
tropy) gas initially at the center to larger radii (see also Li k\ 
|Bryan|[2014| . Second, the dynamical evolution of the AGN 
outflow generates regions of sustained compression (where 
V • r’ < 0), which gives rise to high contrast density inhomo¬ 
geneities. 

As in the model GG, a massive cold filament is created 
along the axis of the jet, with mass ~ 10^ M©, for models 
IG2000 and IG 4 OOO. No gas cools instead for r > 500 pc in 
the case of the more powerful model IG8000 (middle-right 
panel of Figure]^. Evidently, the AGN heating is strong 
enough to completely suppress off-center cooling, at least 
for the timespan probed by our simulation. 


^ Just before the time of the first off-center cooling episode, the 
ratio Gool/Gcc near the cooling region has values in the range 
[20-50]. This is due to the advection of low entropy material from 
the center by the AGN outflow. 
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Figure 7. Temperature maps for the isolated galaxy “single 
event” and i;jet = 2000 km s“^. Top panel: t = 12 Myr; mid¬ 
dle panel: t = 16 Myr; bottom panel: t = 20 Myr. The z-axis is 
horizontal, the i?-axis is vertical. Units are kpc. Notice that the 
color scale differs for every panel. 


We notice that the average total cooling rate in the 
50 Myr following the AGN outflow, for models IGIOOO, 
IG2000, IG4000 and IG8000, is 0.5,0.45,0.32,0.24 M© yr"^ 
(see bottom-right panel of Figure]^- model IGIOOO is not 
shown). These values reflect the (temporarily) suppression 
of cooling due to AGN activity, being the cooling rate of the 
pure cooling flow model Mcooi ~ 0.8 M© yr“^ . 


4.3 AGN cooling in the Low Mass Galaxy 

We describe only very briefly the single outflow models for 
the low mass galaxy. Because of the lower hot ISM density, 
widespread cooling is less efficient for this system. However, 
we shall see below (Section |5.3| that when repeated AGN 
events are considerated, even the LM model undergoes sig¬ 
nificant off-center cooling. 

For the single outflow experiments, we find that only 
LM2000 shows extended cooling, for t ~ 10 — 25 Myr, when 
r\j 10®-3 X 10^ M© are present in the off-center region. Few 
small blobs develop atz~4, i?~0.5 kpc. As usual, a cold 


filament also forms along the z—axis, with mass ~ 3 x 10^ 
M©. 


5 SPATIALLY EXTENDED COOLING: AGN 

FEEDBACK 

The models presented in Section showed that the dy¬ 
namical interaction of a single AGN-driven outflow with a 
smooth and regular ISM generates density (entropy) per¬ 
turbations which lead to (sometime recurrent) widespread 
gas cooling. Understanding these idealized calculations has 
been preparatory to the study of the more realistic problem 
of repeated AGN outbursts, the subject of this Section. 

The intuitive expectation is that in an ISM continuosly 
stirred by the AGN feedback activity the generation of non¬ 
linear perturbations is facilitated, with the result of more 
intense cooling episodes. This is indeed what the simulations 
presented below show. 


5.1 Feedback in Central Galaxies 


In Figure we show the off-center temperature maps for the 
feedback models of the GG galaxy, with outflow velocities of 
2000, 4000 and 8000 km s~^. The dark spots in every panel 
indicate gas cooled to T < 10^ K. In Figure the mass of 
cold gas in the various regions, as defined in the previous sec¬ 
tions, is shown for most models. The first result revealed by 
these images is that a persistent cold ISM phase is naturally 
present, regardless on the feedback details. Of particular in¬ 
terest is model GG-FB16000 — with the fastest outflows 
among our simulations — shown in Figure with the ma¬ 
genta dot-dashed line. In this system the AGN heating is 
strong enough to suppress almost completely the cooling at 
the center of the galaxy (where all the cooling in the classic 
cooling flow picture would occur). Despite that, a strong off- 
center cooling episode happens at t ~ 24 Myr and generates 
a total cold gas mass Mbiob ~ 1.5 x 10® M©. 

The network of blobs and filaments occupies a region 
typically few kpc in size and has mass Mbiob ~ 10® — 10^ 
M©. We are not in the position to calculate the fraction 
of this gas in molecular, neutral or photoionized state. It 
is likely, however, that all these components are simultane¬ 
ously present (see our Discussion below), as supported by 
observations (e.g. Werner et al.|[^14[ [David et al. 2014). 
The spatial distribution of the cold gas generally correlates 
with the direction of the outflows and with the disturbances 
induced in the hot gas. It is worth noting that the bright¬ 
ness of warm gas clumps is expected to decrease with the 
distance from the center, in pace with the decrease of the 
pressure. The aspect and amount of the AGN-triggered cold 
phase are in broad agreement with warm/cold gas observa¬ 
tions in massive ellipticals/galaxy groups ( jGaon et al.|2000 


Werner et al. poU) ). We believe that this is the most robust 


and relevant result of our simulations. The velocity structure 
of the spatially extended cold phase is discussed in Section 

o 

The filament of cold gas forming along the z—axis is al¬ 
ways present in the GG simulations, with mass which varies 
between few 10® to few 10^ M©. The filament has often a 
clumpy appearance, with density variations within it of 2 
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orders of magnitude. The velocity of the cold gas is gener¬ 
ally negative (especially at late times) with typical values of 
100 — 200 km s“^, but outflowing portions of the filament 


are also present. 

The physics of a multiphase ISM and the interaction 
between the phases is very complex, a process poorly un¬ 
derstood even in our own Galaxy (e.g. Cox||200'5 ), and we 
have no chance to give a consistent description here. As for 
the single outflow case, the cold material in our feedback 
simulations condenses in very small structures, usually few 
numerical zones in size. Resolution tests indicate that this 
is true even when the cell size is reduced by a factor of 5 
(that is, adopting AR = Az = 1 pc). The volume filling 
factor in the central region where cold gas is found is always 
very small, usually less than 2-3%. Off-center cold gas usu¬ 
ally spans temperatures in the range 10^ — 10^ K, with only 
a few percent undergoing cooling below 100 K. 

Initially, the cold blobs share the same turbulent ve¬ 
locity field of the hot gas out of which they form. Then, 
they start to fall toward the center of the system, which is 
reached in about a dynamical time (~ 10^ yr); this is the 
typical lifetime of an off-center cloud or filament. 


5.2 Feedback in the Isolated Galaxy 


The isolated galaxy, although more optically luminous than 
the CG, has a much lower amount of hot ISM. As discussed 
earlier on, this leads to a longer tcooi and a larger tcooi/Gcc- 
Correspondingly, the entropy is higher. However, this only 
slightly hinders the spatially cooling process (see also Sec¬ 
tion 4.2). We find that both models IG-FBIOOO and IG- 


FB2000 generate an extended cold ISM phase, whose phys¬ 
ical properties vary in time (see Figures and 10). Total 
masses of cold gas Mcoid ~ 10^ — 10® M© are common in the 
off-center region {R > 500 pc, z > 150 pc). For model IG- 
FB 4 OOO only little gas cools in this volume, Mcoid ~ 3.8 x 10® 
M©, and there is only one episode of extended cooling at 
t ~ 30 Myr. Of course it is possible, perhaps likely, that 
more gas would condense at times > 50 Myr. We did not 
calculate a feedback model for the IG with I’jet = 8000 km 
s“^, but it seems plausible that recurrent outflows with such 
a velocity prevent off-center cooling completely, as the single 


outflow does for the IG8000 model (see Section 4.2) 


The cold gas mass in the z—axis filament is comparable 
with that in the off-center blobs, varying between zero and 
- 7 X 10® M© in models IG-FBIOOO and IG-FB2000, while 
for IG-FB 4 OOO it reaches ~ 7 x 10^ M© at its maximum (cfr. 
top and bottom panels in Figure |^. 


5.3 Feedback in the Low Mass Galaxy 

Results for this model, which has an optical luminosity close 
to the L* parameter of the luminosity function, are relevant 
for a large number of normal elliptical galaxies. Before the 
onset of the AGN, the (hot) ISM mass within the effective 
radius is ~ 1.2 x 10^ M©; within r ~ 2 kpc the ISM is 
slowly inflowing as in the classical cooling flow theory, while 
at larger radii the SNIa heating generates a highly subsonic 
outflow, with average velocity ~ 10 km s“^. Is interesting 
that while gas cooling occurs also where tcooi/Gcc ^ 10, it 
is however limited to the ISM inflowing region (r < 2 kpc. 


Figure [lo|, within which the initial ISM mass is ~ 6.9 x 10® 
M©. We also calculated few models with a higher SNIa rate 
(and therefore a smaller or no inflowing region, although 
the outflow velocity is still very subsonic) to test the rea¬ 
sonable hypothesis that global inflow (except for the outflow 
region, which occupy a small volume of the system) is a nec¬ 
essary condition for widespread cooling. Glearly, non linear 
density perturbations moving outward, toward lower pres¬ 
sure regions, would expand, lengthening tcooi and delaying 
or frustrating the cooling process. In these “windy” galaxies 
we find that no gas cool^ 

In Figure right column, the time evolution of cold 
gas mass in the off-center region and in the r > 500 pc 
region is plotted. For models LM-FBIOOO and LM-FB2000 
an extended cold ISM phase naturally forms within r < I 
(r < 2) kpc for LM-FBIOOO {LM-FB2000), again in an elon¬ 
gated region around the axis. For the high outflow veloc¬ 
ity model LM-FB 4 OOO no gas is able to cool. In this case the 
AGN heating offsets radiative cooling for the (short) time 
covered by our simulation (50 Myr). It is entirely plausible 
that some cooling would happen at later times (see jGaspari 


et al. 2012a), but we did not explore this possibility here. 
Gold gas masses of I — 5 x 10® M© are typical, with periods 
of much lower or no cold gas found in the galaxy. 

As a numerical convergence test we run a model iden¬ 
tical to LM-FBIOOO but with a grid resolution of 2 pc. The 
results are illustrated in the right panels of Figure 0by 
the orange lines. While the detailed pattern of the cold gas 
mass evolution differs, the important result — the genera¬ 
tion of an extended, ~ 10® M© cold ISM phase, even where 
Gooi/Gcc > 10 — is validated. 


6 DUST IN THE HOT AND COLD GAS 
6.1 Central Galaxy 


As mentioned in the Introduction, the dust content of the 
cold gas is a crucial test for its origin. The necessary first 
step is to understand the grain evolution in the hot ISM. 


This has been investigated in detail by [Tsai & Mathews 

(1995 1996 

), Temi et al. (2003), 

Dwek et al.| (|I990|), |Don- 

ahue & Voit| (|I993), Hensley et 

al.| (|20I4|). Stellar ejecta 


continuosly supply dust to the ISM, at a r ate of 
I0“®(M*/I0^^ M©), assuming (5* = I/I50. In Tsai & Math- 


|ews| ( I995| ) and Temi et al. ( 2003| it is argued that stel¬ 
lar ejecta are quickly he ated to the hot ISM temperature 
(see also Mathews 1990) in a timescale < 10® yr (a frac¬ 
tion < 20% of it may remain cold, according to the cal¬ 
culations by Harriott & Bregman 2008| ). Dust grains, as 
they come in contact with ~ lO'^ K gas, are sputtered away 
with a timescale Tgputt ~ 0.03/np Myr (which translates in 
Tsputt < 10® yr in the galaxy core and Tgputt ~ 3 x 10® yr 
at r = 10 kpc), short enough to make grains be disrupted 
essentially where they originate. 

Our model for the dust evolution in the hot gas agrees 
with the aforementioned calculations very well. Dust source 


® In these systems hot gas could in principle still cool through the 
dust assisted cooling process described in [Mathews h Brighenti] 


( 20Q3b| ). 
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Figure 8. Upper panel: time evolution of the cold gas mass in the off-center region {R > 150 pc, 2 ; > 500 pc) for feedback models CG- 
FB2000, CG-FB 4 OOO, GG-FB8000, GG-FB16000 (left panel), IG-FBIOOO, IG-FB2000, IG-FB 4 OOO (middle panel) and for LM-FBIOOO 
and LM-FB2000 (right panel, notice that the time axis extends to 100 Myr in this latter case). Lower panel: as above for the total cold 
gas mass in the off-center region r > 500 pc. 


and sink terms result, for the CG model, in an almost con¬ 
stant ISM dust-to-gas ratio (5 ~ 3 x 10“®. Before the outflow 
starts the total mass of dust present in the gas in the whole 
system is Mdust ~ 7.4 x 10^ M©, of which only ~ 5 x 10^ 
M© reside in the central region within = 10 kpc. These 
numbers are very similar to those of the steady state model 
for interstellar dust in NGC 5044, presented in |Temi et al.| 
(2007), where it has been showed that the observed fluxes 
at 70 jim. and 170 /xm are a factor 50 — 100 higher than ex¬ 
pected if the only source of cold dust were the stellar winds. 
The excess dust in NGC 5044 i s spatially correlated to the 
warm (T ~ 10^ K) gas (see also Gastaldello et al.]|2009 ). 


Gas cooling out of the hot phase initially shares the 
aforementioned low dust content. However, dust can grow 
in cold gas by accretion of condensable elements onto pre¬ 
existing grains, as outlined in Section We find that dust 
accretion is never significant for the cold blobs found in our 
simualtions. The increase in the dust-to-gas ratio of cold gas 
is modest, leading to an average (5coid ^ 10“^ after 50 Myr, 
close to the dust-to-gas ratio of the hot phase. We see little 
variation in the dust-to-cold gas ratio between models. The 
inefficient dust growth originates from both the relative low 
density of the cooled gas clouds (usually < 50 cm“^) and 
their short lifetime. Off-center cold blobs, in fact, are young 
objects, with life span ~ tdyn ^ 10^ yr, after which they 


settle in the galaxy core or in the cold filament along the 
z— axi^ 

Therefore, a robust prediction of our calculations is that 
internally produced, spatially extended cold gas is very dust- 
poor. Evidently, any spatially extended (and therefore re¬ 
cently formed) dusty cold gas in isolated elliptical galax¬ 
ies must have a different origin, as for example cooling of 


dusty stellar ejecta mixed with local hot gas (Mathews & 
BrighentT||2003b Temi et al.||2007 ) or accretion of external 
gas (e.g. Davis et al.||20lT ). On the contrary, cold gas can 
reside in the z—axis filament long enough, and with density 
large enough, to allow dust accretion be effective. At the end 
of the feedback simulations, at t = 50 Myr, the dust-to-gas 
ratio in the filament reaches values of ~ 5 x 10~^. For the 
same reason, dust growth can also be an important process 
in old cold gas located at the very center of a large fraction 
of ellipticals, known to be dusty (e.g. Tran et al.|2QQl [Mar¬ 


tel et al.||2004| |Lauer et al.||2005[ jSinioes Lopes et^iL |2007| 


Alatalo et al.|20r3). 

We remark that the low dust content for the cooled hot 
gas might also lead to a low H 2 and GO content (with respect 
to the Galactic standards), even if the temperature dropped 


® We neglect here the potential role of hot gas rotation, which 
would generate a large disk of cold gas, several kpc in radius (e.g. 
Brighenti fc Mathews|1997a| ). 
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Figure 9. Left column: temperature maps for model CG-FB2000 at four differeut times: t = 12, 20, 28, 48 Myr, from top to bottom. 
The 2 :-axis is horizoutal, the i?-axis is vertical. Uuits are kpc. Time is measured siuce the switch ou of the AGN outflow (which lasts 2 
Myr). Central column: temperature maps for model CG-FBfOOO, showu at 28, 36, 42 aud 48 Myr. Right column: temperature maps for 
model CG-FB8000 at t = 20, 34, 40 aud 46 Myr. Notice that the color scale differs for every pauel. 
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to very low values T < 100 K (see Cazaux fc Spaans|2009 |. 
Furthermore, dust less cold gas would not be effective in ab¬ 
sorbing stellar radiation and would be hardly detectable in 
dust extinction maps. Some indication of dust-poor molecu¬ 
lar gas in NGC 5044 is provided by the ALMA observation 
presented in David et ah (2014). Mainy molecular associa¬ 
tions have been found the inner ~ 2.5 kpc, for a total molec¬ 
ular mass > 5 X 10^ M©. The absorption map presented in 


David et al. (2014) shows dusty patches that do not corre¬ 


late with the molecular structures. If the dust poor cooled 
gas is able to form H 2 and CO in < 10 ^ yr, molecular clouds 
in ETGs would not provide much absorption or extintion, 
just like observed. 

Only the cold z—axis filament is expected to have a 
sizeable amount of dust but, although the formation of this 
feature is entirely reasonable, we are chary to study in detail 
its dust properties here, for the reasons explained above. 


6.2 Isolated Galaxy 

The hot gas dust content in the standard GF model for the 
IG is 2.7 X 10® M 0 . This amount is more sizeable than in the 
corresponding CG model because of the lower average ISM 
density, which makes the grain sputtering process less effi¬ 
cient (Tgputt cx: n“^), especially in the outer galactic regions. 
Within r = 10 kpc Mdust,hot = 1.4 x 10® M©, while within 
re = 6 kpc we have Mdust,hot = 3.6 x 10^ M©. The dust-to¬ 
gas ratio in the hot ISM is now larger than that of the CG, 
rising from Shot ~ 10in the very center, to S ~ 10at 
the effective radius (6 kpc). 

The higher (5hot leads to a higher average dust-to-gas 
ratio for the off-center cold ISM phase, ^coid = 3 — 10 x 10“®, 
for both the single outflow and the feedback models. Higher 
velocity outflows typically generate higher (5coid because gas 
cools at larger radii, where Shot is also larger. The situation is 
therefore similar to that depicted for the GG model: spatially 
extended cold gas originated by cooling of the hot phase 
is predicted to exist in real ellipticals with AGN feedback. 
However, it is expected to be almost dustless; any dust-rich 
cloud or filament has likely been accreted or ejected from the 
nuclear region, where dusty cold gas is known to be present. 

In the 2 :—axis filament, the dust-to-gas ratio never ex¬ 
ceeds 8 X 10“®. The lower value with respect to the CG 
model is mostly due to the lower density in the IG hlament, 
which makes the grain growth process relatively inefficient. 


6.3 Low Mass Galaxy 

The LM galaxy has the lowest average hot ISM density 
among our models. This makes the sputtering time longer 
and the dust content more substantial than in the CG and 
IG systems, with an average dust-to-gas ratio 3.8 x 10“^ 
within 10 kpc and 1.4 x 10“® within at the time of 
the AGN activation. The dust mass in the hot gas is 
Mdust.hot = 1.7 X loyor 1.3 x 10®) M© for r ^ 3.16 (or 
10) kpc, respectively. The cold clumps and hlament gener¬ 
ated by the AGN outbursts are both characterized by aver¬ 
age (5biob,fii ~ 10 ~^ — 10 “^, with the largest values found for 
model LM-FB2000. 
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Figure 11. Left panels: time evolution of cold gas velocity dis¬ 
persion for gas in the off-center region dehned by 0.5 < 2 ; < 5 
kpc, 0.15 < < 3 kpc. The CG-FB, IG-FB and LM-FB 

models are displayed in the top, middle and bottom panels, re¬ 
spectively. For each system we show the feedback model with 
Ujet = 1000,2000,4000 km s“^ with long-dashed green, solid 
black and dashed blue lines (for the CG model we also show the 
r’jet = 8000 km s“^ model with a dot-dashed red line), respec¬ 
tively. Rigt panels: same for the hot gas. 


7 GAS VELOCITY DISPERSION AND COLD 
INFLOW/OUTFLOW 

7.1 Central Galaxy 


It is interesting to investigate the gas velocity dispersion gen¬ 
erated by the AGN outflow. The chaotic motion triggered 
by the outflow expansion and the following buoyant motion 
of low-density regions cause the density perturbations lead¬ 
ing to widespread cooling. Also, it promotes the circulation 
of the metals injected by the SNIa occurring in the galaxy 


(^e.g. Rebusco et al.|2006 Gaspari et al.|2011a|b 2012a). Re¬ 
cent X-ray observations provide estimates on the turbulent 
velocity in clusters and groups, albeit with large uncertain¬ 


ties ( Werner et al.|[2QQ9 Sanders et al.||2011 de Plaa et al. 


2012 ). For NGG 5044 and NGG 5813 the latter authors hnd 
320 < Uturb < 720 km s“^ and 140 < t’turb < 540 km s“^ 
respectively, measured in the central ~ 10 kpc region. 

In the pure cooling flow model, before the onset of the 
AGN feedback, the hot gas 3D, mass-weighted velocity dis¬ 
persion is extremely low, ahot 14 km s~^ (here and in the 
following, the velocity dispersion is measured in the region 
0.5 < z < 5 kpc, 0.15 < R < 3 kpc, corresponding to the 
off-center region where the numerical resolution is higher 


We estimate the mass-averaged velocity dispersion as 
(7^ = where the z—component is cr^ = 

f (vz — Vz)‘^pdVl f pdV and the integrals are calculated on the 
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This value for (Jhot essentially reflects the (smooth) radial 
velocity gradient typical for cooling flows in the region con¬ 
sidered. The single outflow event increases the hot gas tur¬ 
bulence, with typical velocities ahot ~ 15, 35, 60, 100 km 
s“Uor CGI 000 , CG 2000 , GG 4 OOO and CG8000 respectively, 
measured few Myr after the AGN outburst. This chaotic mo¬ 
tion slowly decays, reaching values crhot ~ 15, 30, 30, 40 km 
s“^ at t = 50 Myr. The cold gas, when present, usually shows 
larger velocities, with typical values (Jcoid ~ 80 — 150 km s“^, 
essentially independent on the outflow velocity. These num¬ 
bers compare well with observed cold gas velocity dispersion 
in massive ellipticals ( Caon et al.|2000' Werner et al.|2013 ). 

In Figure we show the velocity dispersion for the 
cold and hot ISM for the feedback models. Again, crhot gen¬ 
erally increases with ujet (although models CG-FB 4 OOO and 
CG-FB8000 are almost undistinguishable). Typical crhot val¬ 
ues range between 50 and 100 km s“^. The obvious differ¬ 
ence with respect to the single outflow models is that crhot 
stays approximately constant in time, as expected for a self¬ 
regulating feedback scenario. 

We also find that crhot is strongly dependent on the po¬ 
sition. As expected, the turbulence is much stronger near 
the jet axis. If we measure crhot in the sector R ^ z (that 
is, far from the jet axis) we get values a factor 1.5 — 2 lower 
than in the sector R ^ z. This implies that any mixing and 
turbulent heating should be much less effective in the region 
near the plane of symmetry z = 0. 

The values for crhot we find are significantly lower than 
those observed. Although estimates for crhot in real galax¬ 
ies are still sparse, this discrepancy may indicate that AGN 
feedback is more powerful than the one considered in the 
present paper (see, e.g. Gaspari et al.|201^ ). Stronger out¬ 
bursts are indeed needed to prevent the excessive total cool¬ 
ing rate observed in our simulations. A further possibility 
to explain the observed large turbulent velocity is the effect 
of merging and accretion into the system. For instance, in 
NGG 5044 a pair of cold fronts are visible within the in¬ 
ner ~ 70 kpc region, likely generated by the accretion of 
a galaxy or a small group (Gastaldello et al. 2009| David 


et al.|[2009 ). The large (~ 150 km s ) peculiar velocity of 


the central galaxy with respect to the mean group velocity 
supports this scenario. 

Also for the feedback models we find that < 7 cold ^ crhot- 
This is in part explained by the tendency of the gas to cool 
in regions of localized stronger turbulence. In addition, once 
cooled to low temperatures, the blobs decouple dynamically 
from the hot gas (although their motion is never totally 
ballisticpi , and this results in a further increase of cTcoid- 
Overall, the recently formed cold gas moves in a somewhat 
chaotic way, with velocities similar to those found by |Gaon| 


et al. (2000). We note here that in our simulations no cold 


gas is dragged to large radii from the central region. This 
indeed must be expected, given the high density and small 
cross section of the cold clouds, which make them almost 
immune to drag by hot gas or entrainment. In models with 


volume defined in the text. The /^—component is computed in a 
similar way. 

The fact that single blobs are not properly resolved in our 
simulations and the effect of numerical viscosity may play a role 
in maintaining a slight coupling between hot and cold gas. We 
did not investigate this in detail. 



-500 0 500 1000 -500 0 500 1000 -500 0 500 1000 


(km/s) 


Figure 12. Distribution of the cold gas in the z-axis velocity 
space for various feedback models at several times. Positive ve¬ 
locity indicate outflowing gas. Every row represents a model, in¬ 
dicated on the left panel of each row (for example in the top row 
is shown the CG-FB2000 run). 


large ujet we do see cold blobs moving outward together with 
the hot gas, but this is because the cold gas forms from the 
outflowing hot gas and initially shares the same motion. 

At any time cold gas clouds with both positive (out¬ 
flowing) or negative (inflowing) radial velocity are present. 
Often gas cools and condenses out from outflowing hot gas, 
resulting in a cold outflow. We note that molecular outflows 


are now commonly detected in AGN hosts (e.g. Gicone et al. 
2014 Tadhunter et al.|2014 ), albeit in systems different than 


the group centered massive elliptical simulated here. In our 
scenario (see also Li & Bryan 2014 Sharma et al. 2012) 
the cold gas is not directly accelerated by the AGN, but 
merely forms in the already forward moving hot gaj^ As 
for the mass and spatial disposition, the distribution of the 
off-center (0.5 < z < 5 kpc, 0.15 < R < 3 kpc) cold gas mass 


In a very recent paper [Costa et a'L| ( |2014| have shown that hot 
gas cooling can generate cold outflows in QSO hosts. 
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in the velocity space, shown in Figurefor several models, 
also varies in an unpredictable way. For most models/times 
there is a slight tendency for the cold gas to prefer outward 
motion, but the effect is subtle. Typical velocities of ~ few 
km s“^, both positive and negative, along the outflow sym¬ 
metry axis are always present. It is interesting that the cold 
gas velocity distribution does not depend in an obvious way 
on the outflow velocity 'Ujet — the model CG-FB16000 does 
not show a high velocity tail as might be expected. 


7.2 Isolated and Low Mass Galaxies 


For these models, the hot gas velocity dispersion before the 
AGN onset is crhot ~ 5 km s“^, somewhat less than in the 
CG case. 

The single AGN outflow in the IG generates chaotic 
motion with crhot ~ 40, 80, 120, 200 km s“^ for IGIOOO, 
IG2000, IG 4 OOO and IG8000 (few Myr after the feedback 
episode, again measured in the region z < b kpc, R < 3 
kpc). The ISM motion dissipates to reach crhot ~ 20 km s~^ 
at t = 50 Myr for IGIOOO and IG2000, and crhot ~ 60 km 
s“^ for IG 4 OOO and IG8000. As for the GG model, IG also 
shows a large velocity dispersion for the off-center cold gas 
(run IG2000)^ with ctcoM ~ 100 km . Again, the main 
reason for this discrepancy is the different spatial regions 
probed by the the two gas phases. While for the hot gas the 
whole off-center volume is used to compute crturb, the cold 
gas is present only in the smaller region {R, z) < (0.8, 2.5) 
kpc, where indeed the hot gas turbulence is stronger and 
comparable to that of the cold gas (~ 100 km s“^). 

The IG feedback simulations exhibit similar charachter- 
istics for the AGN-induced random motion as for the GG-FB 


models (Figure 11). The hot phase 3D velocity dispersion 
varies from ~ 50 to ~ 100 km s“^, while the cold gas shows 
values about twice as those of the hot gas. 

For the low mass galaxy LM, the chaotic ISM motion 
resulting from the AGN activity is similar to that of the IG- 
FB models, and the 3D gas dispersion velocity is shown in 
the bottom row of Figure pT] We note that for the off-center 
cold gas of these LM-FB models ctcoW is somewhat lower 
with respect to the GG-FB and IG-FB runs, a consequence 
of the lower free fall velocity of the LM galaxy. 

For these isolated galaxies the velocity distribution of 
the cold gas (Figure p^ is narrower than for the GG model. 
For most of the time cold outflows are limited to ~ 200 km 
s“^, while the inflowing portion of the cold ISM commonly 
reaches ~ 400 km s“^. 


8 ENERGETICS 


In this work we decided to be conservative regarding the 
strength of the AGN feedback considered. We deliberately 
chose to adopt low power outflows because they are expected 
to be the most common in normal elliptical galaxies. For in¬ 
stance, even in the spectacular archetype of AGN-ISM inter¬ 
action NGG 5044, the pow er associated to the cavity system 
is only < 6 x 10^^ erg s“^ ( David et al.||2009 ). More typical 


cavity powers for m assive ellipticals or groups range from 
10^^ to 10^^ erg s“^ ( Panagoulia et al.|2014| f^ 


We do not show here the detailed evolution of the out¬ 
flow energetics, but limit ourself to a discussion of the basic 
properties of our feedback simulations (see Gaspari et al.| 


2012a for a detailed analysis of more sophisticated models). 


In agreement with the forementioned paper we find that 
the total energy injected by the AGN slightly increases with 
r'jet (although this trend reverts for GG-FB8000). For mod¬ 
els GG-FB1000, GG-FB2000, GG-FB 4 OOO, GG-FB8000 
GG-FB16000 the total outflow energy delivered to the 
ISM at t = 50 Myr is 2.5 x 10^^ 7.0 x 10^^ 1.8 x 10^® 
7.2 X 10^^ and 6.0 x 10^® erg, respectively. As found in 


Gas¬ 


pari et al. ( |2012a ) we also note a definite trend for the duty 


cycle (defined as the ratio between the time in which the 
outflow is active to the total time of the simulation) to de¬ 
crease with increasing r'jet: 2.3 x 10“^, 7.1 x 10“^, 1.9 x 10“^, 
1.5 X 10“^ and 6.7 x 10“^ for the models listed above. This 
is a signature of a self-regulating mechanism. In model GG- 
FB 16000 only a dozen AGN outbursts occur. The aver¬ 
age outflow power, Pjet = F;|etV(60 Myr x duty cycle) is 

6.9 X 10^\ 6.2 X 10^^ 6.0 x 10^^ 3.0 x 10^^ and 6.7 x 10^^ 
erg s~^ for the five considered models. Single AGN outbursts 
show a typical power variation of ~ 5 with respect to the 
mean value. Model GG-FB16000 clearly shows a different 
feedback regime with respect to the lower r'jet calculations. 
The AGN heating operates through a sequence of relatively 
infrequent, very powerful outbursts, separated by quescence 
periods lasting several Myr. While this model is successful 
in reducing the total cooling rate down to acceptable val¬ 
ues, the central ISM often shows temperatures > 2 keV, 
higher by a factor of ~ 2 than those in the center of real 
galaxies/groups. The difficulty for the AGN heating to si¬ 
multaneously halt gas cooling and preserving low central 
temperatures is well known and will not be discussed here 
(see, e.g. Brighenti Mathews|2002 Gaspari et al.||2Q12a ). 

For the isolated galaxy IG the outflows energetics is 
somewhat different with respect of that of the GG models. 
For instance, the total energy injected in the ISM is now 
5.4 X 10^^ 1.1 X 10^® erg for model IG-FBIOOO and IG- 
FB2000, while the duty cycle is 3.3 x 10“^, 8.1 x 10“^ and 

5.9 X 10“^ for IG-FBIOOO, IG-FB2000 and IG-FB 4 OOO. 


9 DISCUSSION 

9.1 The Cold ISM of Elliptical galaxies 

Despite of the accuracy of the simulations presented here, 
we are well aware that many results are subject to quantita¬ 
tive uncertainties. The physics of the AGN feedback process 
still eludes our understanding and has been represented only 
schematically. A source of error which is difficult to evaluate 
in the simulation outcomes originates from the sometimes 
neglected numerical overcooling, together with the intrin¬ 
sic complexity of the cold-hot gas interaction (which is im¬ 
possible to reproduce in current galaxy-wide calculations). 


It should be noted here that these estimates for low mass sys¬ 
tems might be inaccurate, if very deep Chandra observations are 
not available. Furthermore, the fraction of outflow power dissi¬ 
pated in shocks is poorly known. 
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Notwithstanding this cautionary note, we believe that a solid 
result of our study (see also Brighenti, Mathews and Temi, 
2015) is that AGN activity triggers spatially extended gas 
cooling from the hot ISM phase. This happens for a wide 
range of outflows parameters, a fact which makes this pre¬ 
diction robust despite the uncertainties in the AGN feedback 
procesj^ Of course, the net effect of the energy injected by 
the AGN is to reduce (or suppress) the total cooling rate 
(e.g. Gaspari et al.|2012a ). 


The off-center cold gas has an erratic distribution, fol¬ 
lowing the cycle of the AGN feedback, activated (by design) 
by the random presence of cold (dropped-out) gas near the 
black hole (Pizzolato & Soker 2005 2010| Gaspari et ah 


2013 2014). The exact amount of cold gas depends on the 


outflows parameters (velocity, especially) and is therefore 
somewhat uncertain. In our CG model (tailored on NGG 
5044) we find that ~ 10^ — 10® M© of off-center cold gas 
are typically present, depending on the feedback parame¬ 
ters. For the isolated galaxies IG and LM the amount of 
the cold gas is reduced, but still significant: ~ lO'^ - 10® 


M©. These values compare well with observations (Mac- 


chetto et al.|1996] Gaon et al.|2000 David et al.|2014 ). Espe¬ 
cially for the isolated systems, there are times when the cold 
gas is absent or negligible (see Figure]^. Indeed, that the 
mass budget of the cold ISM phase undergoes unpredictable 
variations (cfr. Figure]^, again in pace with the AGN activ¬ 
ity. It seems therefore plausible that AGN-induced cooling 
can explain, at least in part, the large fraction of intermedi¬ 
ate/massive ellipticals with detectable Ha emission or even 
molecular gas, if H 2 and CO can form fast enough in the 
dust-poor cooled gas. 

Our simulations cannot follow the evolution of the 
cooled gas with the accuracy necessary to calculate what 
fraction of it ends up in molecular form. H 2 can form by 
both grain surface and gas phase reactions, with the dust 
grain processes dominant fo r gas metallicity larger than 


Z = I0 — 10^ solar (e.g. Gazaux & Spaans 


2009 


and 


references therein). The cold gas in our simulations has ~ 
solar metallicity, but the dust-to-gas ratio is (5 ~ 10“^. This 
latter value is similar to the dust -to-gas ratio of the mode l 


with metallicity Z = 10 solar in Gazaux & Spaans (2009), 


which we can therefore adopt to make order-of-magnitude 
considerations. Using the H 2 formation rates given in their 
Fig. 6, we find that the typical H 2 formation timescale is ~ I 
(50,2000) Myr for gas density ~ 10® (10^, 10) cm“®. Unfor¬ 
tunately, for the highest quoted densities, these numbers are 
close to the lifetime of our simulated cold clouds, therefore 
no strong conclusions can be drawn about the presence of 
molecular gas in the blobs. For n = 10 cm“® no H 2 is ex¬ 
pected to form. It is unclear which density the gas cooled 
to T < 100 K would reach — our simulations do not have 
the resolution to properly resolve the structure of the cold 
clouds. Moreover it seems likely that magnetic or turbulent 
support would be important in regulating the density pro¬ 
files of the clouds. Thus, in order to properly investigate this 
topic, a better calculation for the thermal and dynamical 
evolution for the clouds is needed, coupled with a H 2 (and 


Numerical experiments with different feedback methods, such 
as jets and cavities generated by cosmic rays (Perrotta et al., in 
preparation), yield qualitatively similar results. 


CO) formation model which takes into account the chemical 
properties of the cooled ISM and the appropriate dissociat¬ 
ing UV field. This is beyond the scope of the present paper, 
although such a calculation would be crucial to properly in¬ 
terpret the molecular associations detected in the central 
galaxy NGG 5044 by David et al. (2014). 


Davis et al. (2011) found that most galaxies show kine¬ 


matically misaligned cold gas. This is usually interpreted as 
evidence that gas has been accreted, through (minor or ma¬ 
jor) merging or from smooth cosmological accretion. That 
is, that cold gas has an external origin. Another possibil¬ 
ity is the scenario proposed by Lagos et al. (2014), where 
misaligned cold gas is generated internally by cooling of a 
twisted hot gas halo. Undoubtedly these processes occur and 
may well explain the bulk of the misaligned cases. However, 
it is interesting to note that AGN cooling also predicts kine¬ 
matically decoupled cold gas-stellar systems. In fact, cold 
gas forms and moves preferentially along the outflow axis, 
which may well be the rotation axis of the stellar body of the 
galaxy. Therefore, looking to the system from a generic view¬ 
ing angle, stars and cold gas would indeed show a kinematic 
misalignment, even though the gas origin is internal. Of 
course, the cold gas in this simple scenario would likely show 
irregular appearance and kinematics — smoothly rotating 
misaligned gas can be internally generated only through the 
model described in Lagos et al. (2014). 

Thus, AGN can play a double role in the galaxy life. It 
strongly decreases the gas cooling process, solving the cool¬ 
ing flow problem and quenching the star formation, causing 
the migration of the galaxy to the red sequence. However, 
the dynamical interaction with AGN outflows keeps the ISM 
“alive”, giving rise to the episodic presence of a significant 
cold phase and possibly triggering sporadic, weak star for¬ 
mation events which could explain in part the recent SF 


activity commonly detected in ETGs (Trager et al. 2000 
Kaviraj et al.|2007 ). While minor mergers are an attractive 


explanation ( Kaviraj et al.|[2009 ), we point out that AGN- 
induced hot gas cooling is a natural, internal process which 
also potentially leads to low-level SF. 


9.2 Ha emission from warm clouds 

Gold blobs, although they share similar properties, are made 
up of gas which has cooled down to different temperatures 
(1.8 < /o^'Tbiob ^ 4). Here we provide a rough estimate of 
the expected Ha luminosity of the spatially extended cold 
gas, which in turn depends on the amount of the ionized 
cold ISM. Golumn density of blobs, calculated as N = {n)l, 
where (n) is an average of the blob number density and I is 
its linear size, spans the range 3 x 10^^ ^ ^ ^ 3 x 10^^ cm“^. 
Their number density is in the range 0.4 < n < 50 cm“® and 
their (smallest) characteristic size is usually no more than 
few tens of parsec. This column density is usually larger 
than the maximum column density iVion that can be pho- 
toionized by the local ionizing photon field, which includes 
the contributions of both the stellar radiation (usually the 
dominant one at the clouds location) and the intergalactic 
background. 

The mass fraction Mion of ionized gas in blobs has been 
estimated in order to derive the actual contribution of cooled 
gas emitting in Ha band. Assuming that blobs are spherical, 
the problem for a single blob reduces to an inverse Stromgren 
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sphere: AttRs^ /S.RiF a — AttRs^ riryC^ 
of the blob (or the smallest size of an elongated clump), Ai? 
is the thickness of the actually ionized shell, n and the 
density of the gas inside the cool clump and the local ion¬ 
izing photon density respectively, c the speed of light. Here, 
the number of recombinations expected in the spherical shell 
has been equated to the number of ambient incident photons 
which photoionize the outer layer of the blo b; the recombi - 
nation coefficient is a ~ 4 x 10“^^ cm^ s“^ ( Spitzer||l978 ). 
The ionizing photon density n^^r) has been estimated by 
following Brighenti &; Mathews| ( |1997b ), who considered an 
elliptical galaxy of stellar mass 7.5 x 10^° Mq and effective 
radius 1.72 kpc, by scaling their n^{r) for the different stel¬ 
lar masses and effective radii of our galaxies. The resulting 
ionizing photon density profile writes cm“^, 

where the constant A = 10“^, 3 x 10“^, 2.4 x 10“^ for CG, 
IG and LM respectively. 

Solving for Ai^, the thickness of the actually ionized 
layer of a typical blob which is located at a certain dis¬ 
tance from the center at a given time and which has a 
certain density is estimated. From the value of AR ^ 
(20, 70, 58) pc (for models GG, IG, LM) it is pos¬ 

sible to calculate the ionized mass of the blob: 

Ai7 

Mion ~ AttRs A_R Pbloh ^ ‘AMblob 5 (0) 

JLs 

where puob is the density corresponding to the aforemen¬ 
tioned n and Rs — 30 pc. Unfortunately, we find that the 
fraction of ionized gas in a single clump varies widely, de¬ 
pending on its density and location in the galaxy, spanning 
the range 0.2% Muob ^ Mion ^ 50% Muob^ 

It is difficult and likely not worthy, given the uncertain¬ 
ties, to calculate the percentage of ionized mass of every sin¬ 
gle blob. However, density and temperature of a large num¬ 
ber of examined clumps span the aforementioned ranges; the 
simulation time does not play a crucial role in determining 
the density or the temperature of appearing (just cooled) 
blobs; on the contrary, the evolution of blobs, when it can 
be followed, shows that old blobs are denser than short¬ 
lived ones and that they are often the result of the merging 
of ligther clumps. 

A further difficulty is that the simulated blobs are 
under-pressurized with respect to the ambient gas, some¬ 


times by a factor of 10 — 100 (see also Li & Bryan 2014). 
If temperatures of blobs are thought to be reliable, the ex¬ 
pected pressure equilibrium would produce a correspond¬ 
ing increase in the density (neglecting magnetic support 
This, in turn, leads to a smaller linear size of the blob. The 
combined effect of the increase of puob and the drop in I low¬ 
ers the mass fraction of ionized gas by a factor ~ 5 x 10~^ 
when a typical blob of given mass is thought to have a two 
order of magnitude larger density. The aforementioned in¬ 
crease in pbiob would lead to the enhancement of dust growth 
in blobs, too, since p growth cx^ P (see section]^. However, a 
cloud is expected to reach pressure equilibrium after a sound 
crossing time Tsc\ Tsc spans the range 2 — 20 Myr when blobs 
of a typical size of 30 pc and of different temperatures are 


Self-gravity, whose effect is not included in the code, cannot 
physically play a crucial role in attempting to compress blobs 
and enhance their density, the smaller size of a typical blob being 
usually far shorter than its Jeans length. 


considered. The long Tsc might be the responsible for the 
lack of pressure equilibrium in coldest blobs, most blobs be¬ 
ing younger than their Tsc and more short living. Observa¬ 
tions indicate low thermal pressure in the ionized gas. For 
instance, in the central regions of M 87 the [S H]A6716/6731 
line ratio reveals electronic number densities (rie — 25 cm“^) 
in the emission line gas (Werner et al. which are lower 

than the ones obtained if pressure equilibrium is supposed 
[ue — 100 cm“^ in gas with T ~ 10^ K if a typical number 
density of ne,hot — 0.1 cm“^ in the That ~ 10^ K surround¬ 
ing medium is adopted). 

These considerations suggest to not attempt a rigorous 
calculation of the Ha luminosity in our simulated galax¬ 
ies. Nevertheless we are able to give interesting upper limits 
in Lhcx generated by the warm non-nuclear gas. We rerun 
several models similar to those described in the previous 
sections, but this time setting T = 10^ K as lower limit in 
the gas temperature and assuming that all the cold gas is 
ionized. In this way we can calculate firm upper limits for 
Lhcx- For the GG-FB system we find for the off-center gas 
an average value LnaMoh ~ 5 x 10^^ erg s“^, with typical 
variations of ~ 5 among the various GG models. The z-axis 
filament can be much brighter, LhcxJU ~ 10^° — 10^^ erg 
s“^, but most of its gas is indeed expected to be neutral or 
molecular. The off-center warm gas in IG and LM models is 
fainter than for the GG, an obvious result of the less sub¬ 
stantial cold ISM phase in these systems. The upper limits 
for the IG-FB (when off-center cold gas is present — see 
Figure]^ vary in the range LHa, blob ~ 10^^ — 10^^ erg s“^, 
while LM-FB never exceeds ~ 5 x 10^^ erg s“^. As before, 
the upper limits for the filament can be one or two order of 
magnitude brighter that those for off-center blobs. 


10 CONCLUSIONS 

In this paper we have presented hydrodynamical models of 
AGN outflows expanding in the hot ISM of elliptical galax¬ 
ies. We have designed these numerical experiments to test 
the idea that AGN feedback, while strongly quenching the 
total cooling rate, can nevertheless stimulate some degree of 
hot gas cooling at a few kpc distance from the galaxy center 
(see also Brighenti fc Mathews|2002 Gaspari et al.|2012a|b 


[Sharma et al.|20T2 McGourt et al.|2012 Gosta et al.|2014 ). 

This process can generate a persistent cold ISM phase in 
early-type galaxies. We explore the feasibility of this sce¬ 
nario in galaxies of different mass and hot gas content, rang¬ 
ing from X-ray bright, group centered ellipticals to galaxies 
of intermediate mass with a relatively thin hot gas halo. 

We have deliberately chosen a simple AGN feedback 
scheme. We assume that AGN injects energy and momen¬ 
tum through non relativistic outflows, activated by the gas 
cooling in the inner 100 pc (see e.g. Gaspari et al.|201^ ). We 
focus here on relatively weak feedback, likely the most com¬ 
mon regime in normal galaxies. The most significant results 
of our investigation can be summarized as follow: 

i) The interaction between AGN outflows and the ISM 
generates non-linear density perturbations in regions of con¬ 
verging flows. These perturbations trigger localized hot gas 
cooling at distances of several kpc from the galaxy cen¬ 
ter. This confirms the early results of [Brighenti Mathews] 
(2002) for the GG model and extend them to elliptical galax- 
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ies of lower mass (IG and LMruns). Li & Bryan (2014) have 
shown that this process occurs also in BCGs at the center 
of massive clusters. 

a) The cooling process and the amount of off-center 
cooled gas are not sensitive on the parameters of the out¬ 
flows. This is a key result, given the uncertainties of the AGN 
feedback mechanism and makes our findings quite robust: 
AGN outbursts seem to stimulate large distance cooling no 
matter what. 

in) Off-center cooling occurs when tcooi/tff ^ 70 (in the 
unperturbed ISM), but we do not expect that this is a strict 
threshold. Instead, the spatially extended cooling process, 
no more closely related to thermal instability when the in¬ 
homogeneities are strongly non-linear, depends on the local 
(perturbation) value of tcooi/^ff, that is on the perturbation 
density contrast (see also |Joung et ah' 2012). Our calcula¬ 
tions suggest that even moderate AGN feedback is sufficient 
to generate density (entropy) inhomogeneities which cool lo¬ 
cally. 

iv) The amount and the spatial distribution of the off- 
center cold ISM is erratic (see Figure]^. For group-centered 
galaxies (model CG), the extended cold phase has a typi¬ 
cal mass of 10® — 10^ M©, with time variation of a factor 
of ~ 5. For isolated galaxies {IG and LM) periods of no 
off-center cold gas give turns to episodes of significant cold 
and extended ISM presence, with mass ~ 10^ — 10® M©. 
A possible spin-off result concerns hot gas halos in massive 
spiral galaxies. Given the similarity between detected hot 
coronae around disk galaxies and the hot ISM in our LM 
model, it appears conceivable that any outburst from the 
central black hole (perhaps similar to the one responsible 
for the Fermi bubble) could indeed trigger cooling of the 
hot corona, forming cold clouds which will accrete on the 
galactic disk. 

v) The cold gas initially shares the same kinematics of 
the hot out of which it forms. Dispersion velocities in the 
range 150 — 250 km s“^ for model GG, 100 — 200 km s“^ 
for IG, 50 — 100 km s“^ for LM are typical. This compares 
well with measurements for the emission line component in 
ETGs (e.g. Gaon et ffi^|2000 ), while is generally lower than 
the values given in Werner et al. (2014) for [GII] emitting 
gas in a sample of massive ellipticals. 

vi) Some gas condenses out from outflowing hot gas, re¬ 
sulting in a cold outflow movin g with characteristi c velocity 
of a few lOOs km s“^ (see also Gosta et al. ]|2014| |. We also 
notice that cold gas moving mainly close and along the out¬ 
flow axis can appear kinematically misaligned with respect 
to the stars. 

vii) We follow the evolution of the dust in both the 
hot and cold gas. Dust is injected in the hot phase by the 
stellar mass loss, it is destroyed by thermal sputtering and 
can grow in the cold clouds by accretion of condensable el¬ 
ements. We find that the off-center cold gas always has a 
very low dust-to-gas ratio, of the order of (5coid ~ I0“® for 
the GG and (5coid ~ I0“® — 10“^ for the lower mass sys¬ 
tems. Therefore, a strong dust deficiency is the distinctive 
sign for cold gas originated internally from hot halo cooling. 
However, we acknolewdge some degree of uncertainty for the 
latter conclusion. If cold gas densities were grossly underes¬ 
timated in our simulations (see Section 9.2) dust could grow 
fast enough to restore a standard dust-to-gas ratio. 

via) It is difficult to predict the fraction of cold gas 


which is ionized by the galactic ionizing radiation. However, 
upper limits on the radiation emitted by the off-center 
cold gas range from ~ 5 x 10^^ to ~ 10^®, with the higher 
luminosity appropriate for the GG system. 

ix) In the Appendix we show that the cold gas masses 
found in our simulations are likely overestimated by a fac¬ 
tor of a few. This is due to numerical mixing and the fol¬ 
lowing overcooling — a phenomenon important when small 
cold clouds move through a much hotter medium, sometimes 
overlooked in the literature. 
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APPENDIX A: NUMERICAL MIXING AND 
OVERCOOLING 

Numerical mixing at the interfaces between cold and hot 
gas can cause large errors in the computed cooling rate — 
not always recognized in the literature. This is especially 
true when the cold gas clouds are resolved in only a few 
zones and when the difference in temperature and density 
between the two phases is sizable. These are exactly the cir¬ 
cumstances of our simulations, where small blobs of T < 10^ 
K gas move through ~ 10^ K hot atmosphere. Numerical av¬ 
eraging at the surfaces separating cold and hot gas results 
in zones of gas at intermediate density and temperature. 
Radiative cooling can be either overestimated and under¬ 
estimated, depending on the flow parameters, with the for¬ 
mer case dominating in our simulations. Brighenti, Math¬ 
ews & Temi (2014) also discuss the magnitude of this errors 
in simulations similar to those presented here. They argue 
that adding a temperature dependent “dropout” term to the 
continuity equation, which removes from the grid the cur¬ 
rently cooling gas, largely reduces the spurious overcooling. 
The trade off of this procedure is, of course, the inability to 
follow the subsequent dynamics of the cooled gas; only its 
amount and the location of its formation can be estimated. 
This term has already been introduced in Section 2.2.2, see 
Brighenti, Mathews & Temi (2014) for more details. The 
term added to the continuity equation reads: 

(§?) = l-iT/Gf] 

(Al) 

and a corresponding term is needed to the energy equa¬ 
tion (see Brighenti, Mathews & Temi 2014). Here we choose 
Tc = 5 X 10^ K. These dropout terms are negligibly for 
T > Tc. By targeting only gas with T < Tc means that we 
remove gas that is currently cooling and will reach T < 10^ 
K in a very short time. That is why this method produces 
accurate cooling (dropping out) rates, providing Tc < 10® 
K (see [Brighenti fc Mathews||2002| ). The precise value of qo 
hardly matters and we choose go = 2 to match our previous 
calculations. 

Needless to say, this technique to alleviate numerical 
overcooling is not perfect. A sensible way to gauge the mag¬ 
nitude of this error is to calculate similar models with and 
without the dropout term (acting in the whole numerical 
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Figure Al. Left column: density maps for the “one-zone” experiment at three different times. In the left panel the “x” indicates the 
original location of the single cold zone. The z-axis is horizontal, the i?-axis is vertical. Units are kpc. It is evident that the small cloud 
induces a “snowball” cooling effect as it moves through the grid. Notice that the color scale differs for every panel. 


grid), and compare the amount of cooled (dropped out) gas 
in the two cases; the physical truth lies in between these two 
extrema, with the dropout value likely to be more accurate. 


Al “One cold zone” experiment 

Before to make this comparison for the CG model, we show 
with a suite of idealized, somewhat pedantic, numerical ex¬ 
periments the occurrence and the effect of the overcooling. 
We follow the evolution of a one-zone cold cloud, located in 
the same gas atmosphere of the CG model at (i?, z) = (7, 7) 
kpc. The temperature of the zone is set to 10^ K and its 
density is calculated assuming pressure equilibrium with the 
surrounding hot gas. The numerical resolutions adopted for 
the calculations are Ai? = Az = 20 or 10 pc. The dropout 
term above is not used in this experiment, and the cold gas 
is simply advected under gravity toward the galaxy center. 
The cloud evolution for the radiative case, resolution 20 pc, 
is shown in Figure [Al] 

At t = 20 Myr the mass of the cloud has grown by a 
factor 2.4, as can be seen in Figure [A^ (solid line), where we 
show the cold gas mass evolution (normalized to the initial 
cold mass). This increase in the cold gas mass is spurious 
and it is due to numerical overcooling. While this mimics 
some physical processes (for example the condensation for 
clouds larger than the Field length, when thermal conduc¬ 
tion is present), we stress that in the hydrodynamical equa¬ 
tions solved no transport terms (which deal with microscopic 
fenomena such as diffusion) are present, therefore any incre¬ 
ment of the cold gas mass must be regarded as spurious. At 
the end of the calculation (50 Myr), the cold gas has grown 
by a factor > 40. 

If the radiative cooling is turned off, numerical mixing 
consumes the cold cloud, in a time t ~ 10 Myr (Figure A2 


pc). The first shows that a (radiative) single zone cloud in¬ 
creases its mass by a factor > 100 (Figure A2 dotted line). 


dashed line). After this time the original cloud gas has been 
mixed with hot gas and raised its temperature to values 
> 5 X 10^ K, our threshold for the definition of cold gas. 
We run two more simulations at higher resolution (10 


In the second, high resolution experiment we consider a 4- 
zones cloud (which is then physically identical to the single¬ 
zone cloud at 20 pc resolution). In this case, illustrated in 
Figure [A^ with a long-dashed line, the cold gas mass growth 
is almost identical to the model described by the solid line. 
This implies that a modest increase in the numerical reso¬ 
lution likely does not solve the numerical mixing problem. 

These few simple calculations, although extreme, are 
relevant for our calculations. It is our experience that regions 
of cold (~ 10^ K) gas only one or a few zones across can arise 
naturally in radiating gas dynamics (Koyama & Inutsuka 


2004). Thus, small elements of cold gas can trigger large 


unintended numerical overcooling due to unphysical mixing 
with much hotter gas. Conversely, in very slowly cooling or 
adiabatic flows, some degree of numerical overheating can 
occur. The physical truth lies between these two extremes. 


A2 Numerical overcooling in the CG-FB models 

In this section we address the effect of numerical overcool¬ 
ing in some CG-FB models. We also performed experiments 
for IG-FB calculations, with similar results. For the sake of 
conciseness, we present results for only two runs of CG-FB 
with jet velocities r'jet = 4000 and 16000 km s“^. These 
calculations are identical to those described in Section [5T] 
but the term in equation |A1| (and an analogue term in the 
energy equation) is activated on the whole numerical grid. 
In Figure |A3| we show the cumulative mass dropped out 
in the off-center region and compare it with the cold gas 
mass found in the same region for the respective standard 


models (see Section 5.1). The masses which are portrayed 
in Figure [A^ have different meaning. For the standard mod¬ 
els, the displayed mass represents the instantaneous cold 
gas mass present in the off-center region: this cold gas then 
flows toward the center, disappearing from that region. For 
the dropout runs it is instead the total gas mass which has 
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time (Myr) 

Figure A2. Time evolution of the cold gas mass for the one- 
zone (and one four-zones) experiments. For every run the mass 
is normalized to the initial cold mass. In radiative simulations 
numerical mixing causes overcooling, while in the adiabatic cal¬ 
culation the spurious mixing erases the cloud in a short time. 


{R,z) ~ (1.5,4) kpc. These blobs, initially only a few tens 
pc in size, grow in size and mass by numerical mixing as 
they move through the grid, exactly as depicted in Figure 
|A1[ That is, overcooling has been triggered. The mass dis¬ 
crepancy can be as large as one order of magnitude, and will 
be not reduced by increasing the numerical resolution by a 
factor of a few. Therefore, the cold gas mass values quoted in 
the previous sections are uncertain by (at least) a factor of a 
few. Although we argue that the dropout masses are closer 
to the physical truth than the cold gas masses calculated in 
the standard fashion, we focused on the latter estimates in 
the text, mostly for the sake of tradition. However, we stress 
the need to evaluate, whenever possible, the uncertainty in 
the cooling gas mass due to numerical overcooling, using a 
drop out procedure or other appropriate experiments. 



Figure A3. Time evolution of the cold gas mass for two CG- 
FB models. The solid blue line indicates the same CG-FB 4 OOO 
model shown in Figure In this calculation no dropout term 
has been used. The blue dotted line represents the same model, 
but with the term IMl adopted. The two magenta lines show the 
standard (dashed line) and the dropout (dot-dashed) version of 
GG-FB16000. 


been cooled and removed from the grid since t = 0. That is 
why it is a monotonically increasing value. 

The effect of the overcooling in the case of no-dropout 
is clear, especially for the i;jet = 16000 km s~^ model. In 
the temperature map, not shown here, it is evident as a 
large number of small cold blobs form at t ~ 25 Myr for 
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